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Abstract —This paper describes a novel framework for computing geodesic paths in shape spaces of spherical surfaces under an 
elastic Riemannian metric. The novelty lies in defining this Riemannian metric directly on the quotient (shape) space, rather than 
inheriting it from pre-shape space, and using it to formulate a path energy that measures only the normal components of velocities 
along the path. In other words, this paper defines and solves for geodesics directly on the shape space and avoids complications 
resulting from the quotient operation. This comprehensive framework is invariant to arbitrary parameterizations of surfaces along 
paths, a phenomenon termed as gauge invariance. Additionally, this paper makes a link between different elastic metrics used in the 
computer science literature on one hand, and the mathematical literature on the other hand, and provides a geometrical interpretation 
of the terms involved. Examples using real and simulated 3D objects are provided to help illustrate the main ideas. 
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1 Introduction 

In this paper we seek a framework for analysing shapes 
of a certain class of 3D objects. Although the general goal 
in shape analysis is to develop tools for full statistical 
analysis - statistical averaging, finding principal modes 
of variations in a population, and shape classification, 
we restrict to more basic goals of quantifying shape 
differences and generating deformations. While there 
have been many efforts in shape analysis of 3D 
objects, the problem is far from solved and the current 
solutions face many technical and practical issues. For 
instance, many general techniques for shape analysis 
rely on quantifying shape differences by spatially 
matching geometric features across objects. Therefore, 
it becomes important to establish a correspondence 
of parts between objects, i.e. which part in one object 
corresponds to which part in the other? This was an 
important bottleneck in a majority of previous efforts 
on 3D shape analysis where the correspondence (or 
registration) of objects was either presumed or solved 
as an independent pre-processing step. More recently, 
there has been progress in establishing frameworks that 
formulate the registration and comparison problems 
jointly. These newer frameworks, using techniques 
from differential geometry, focus on shape analysis of 
parameterized surfaces and treat the problem of shape 
comparison as the problem of computing geodesic paths 
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in shape spaces under a chosen metric. Here shapes 
are compared using a Riemannian metric on a pre-shape 
space T consisting of embeddings or immersions of a 
model manifold (like the sphere, or the disc) into the 
3D Euclidean space Two embeddings correspond 
to the same shape in if and only if they differ by 
an element of a shape-preserving transformation group, 
such as rigid motion, scaling, and reparameterization. 
The shape space is therefore the quotient space of the 
pre-shape space by these shape-preserving groups. 
If the Riemannian metric on the pre-shape space is 
preserved by the action of the shape-preserving group 
then it induces a Riemannian metric on the quotient 
space. The construction of geodesics in shape space 
provide optimal deformations between surfaces and is 
a very important tool in statistical analysis of shapes. 
Interestingly, the problem of registration is handled 
using parameterizations of surfaces such that the points 
denoting the same parameter values on two objects are 
considered registered. 

While these geometric ideas are powerful and com¬ 
prehensive, there are two important issues that one 
needs to deal with: (1) the choice of Riemannian metric 
to define geodesics, geodesic lengths, and the eventual 
shape metric, and (2) the task of computing geodesic 
paths between arbitrary shapes. In terms of the first 
issue, the choice of a metric, an important requirement 
is that the metric should be invariant to action of the 
reparameterization group, to enable a well-defined dis¬ 
tance on the eventual quotient space or the shape space 
of surfaces. There is a related requirement for the shape 
analysis to be invariant to parameterizations of objects 
since parameterizations are only artificial impositions 
designed to help navigate along objects. The physical 
intuition we have is that shape tools, such as the de¬ 
formation (path or geodesic) from one shape to another. 



IEEE PATTERN ANALYSIS AND MACHINE INTELLIGENCE 


2 



Fig.l. Two paths in T with the same sequence of shapes but with different reparameterizations between the corresponding shapes. 


are physical processes that are independent of the way 
surfaces may be parameterized. These dual requirements 
rule out the use of commonly-used quantities such as the 
norm on the space T directly. In terms of the second 
issue, the lack of standard metrics makes it complicated 
to compute geodesic paths even when the underlying 
manifold is a vector space, and one needs numerical 
algorithms for approximating geodesic paths. Next, we 
present a summary of the past work on these two issues 
and outline motivations for the current paper. 

1.1 Motivation and Past Work 

Our goal is to develop tools for analyzing shapes of 
two-dimensional surfaces with certain local constraints 
(smoothness, no-holes, etc). The main difficulty in com¬ 
paring shapes of such surfaces is that there is no pre¬ 
ferred parameterization that can be used for registering 
and comparing features across surfaces. Since the shape 
of a surface is invariant to its parameterization, one 
would like an approach that yields the same result 
irrespective of the parameterization. 

Furthermore, we are not only interested in the com¬ 
parison and matching of two shapes, but also in the 
deformation processes that may transform one shape 
into another, i.e. metamorphosis. To be physically mean¬ 
ingful, the evolution from one shape to another should 
be independent of the way surfaces may be parameter¬ 
ized. Our approach to shape analysis presented in this 
paper was therefore initiated by the following question : 
W/zflf is the natural framework where one can measure de¬ 
formations of shapes independently of the way shapes are 
parameterized? As a motivating example, the sequence 
of shapes displayed in Fig. (bottom) denotes a path 
where a horse is transformed into a jumping cat. During 
the transformation process, only the change of shape, 
drawn in the bottom line as a sequence of blue surfaces, 
is relevant to us. How the surfaces may be parameterized 
during the metamorphosis has no importance in our 
context. To emphasize this idea, two paths of parameter¬ 
ized surfaces corresponding to the same transformation 
process are displayed in the top two rows. We would like 


a framework where the physical quantities measured on 
the path of shapes, such as its length or its energy, are 
independent of the parameterizations of surfaces along 
the transformation process. In particular, in Fig. the 
two paths of parameterized surfaces corresponding to 
the same transformation process should have the same 
length. Note that the surfaces along the second path are 
obtained by applying a different reparameterization at each 
time step to the surfaces along the first path. 

Let us emphasize that we are not only interested 
in how far the horse and the jumping cat are from 
each other, in other words in a quantity like a distance 
measuring the minimal cost needed to deform the horse 
into a cat. But, given a metamorphosis between these 
two shapes, we are also interested in measuring its 
length on one hand, and its energy on the other hand, 
independently of the parameterizations of the transfor¬ 
mation process that may have been used to create this 
metamorphosis. Recall that the length of a path is the 
integral of the norm velocity function with respect to 
time and has the dimension of a distance. The energy is 
the integral of the square of the norm velocity function 
with respect to time, hence has the dimension of the 
square of a distance divided by time. 

Let us now summarize past work on related sub¬ 
jects. The initial set of papers developed algorithms for 
geodesic deformations between surfaces while using the 
given registration of points. They compute geodesics 
between shapes, under isometric deformations, while 
assuming the registration (or parameterization) as given. 
Windheuser et al. ||T| proposed to find a geometrically 
consistent matching of 3D shapes which minimizes an 
elastic deformation energy but use a linear interpolation 
between registered pairs of points in to compute 
geodesic paths. Another paper by Kilian et al. ||2| rep¬ 
resents parameterized surfaces by discrete triangulated 
meshes, assumes a Riemannian metric on the space of 
such meshes, and computes geodesic paths between 
given meshes. The main limitation here is that it assumes 
the correspondence between points across meshes. That 
is, we need to know beforehand which point on one 
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mesh matches with which point on the second mesh. 
The same limitation holds for the paper by Heeren et 
al. m also. In contrast, we would like to remove the 
reparameterization variability so that different surfaces 
with the same shape but different parameterizations 
have zero distance between them. 

Motivated by progress in shape analysis of curves lH, 
0, Kurtek et al. 0,0 introduced a new representation, 
termed a q-map of surfaces such that the distance 
in this representation space is invariant to simultaneous 
reparameterizations of surfaces. For convenience of the 
reader, we recall the definition of the q-map but we 
will not use it in the present paper. Let / : ^ 

denote a smooth parameterized surface and T be the 
set of such surfaces. Then, this g-map is given by / g 
where g(s) = ^yr{s)f{s) and r{s) is the area multipli¬ 
cation factor of / at 5 G §^. They defined a Rieman- 
nian metric on the space of parameterized surfaces by 
pulling back metric under the g-map, and used a 
path-straightening algorithm to compute geodesic paths 
between given surfaces in a pre-shape space. This path¬ 
straightening is an iterative algorithm that updates an 
arbitrary initial path using the gradient of the energy 
function mentioned above, until the path converges to 
a geodesic. The energy gradient is approximated nu¬ 
merically using an (approximate) finite basis for T. To 
remove the effects of original parameterizations, and 
to obtain geodesics in the shape space, they solve for 
an optimal reparameterization of one of the surfaces, 
under the same energy. There are several other papers, 
including |8|, that focus exclusively on the task of finding 
optimal correspondence between 3D objects, either using 
physically-motivated energies or Riemannian metrics. 
Due to the use of gradient-based searches, these methods 
and previously mentioned papers do not guarantee a 
global solution, either for geodesics or for registration. 
In path-staightening, however, it can be shown that a 
path that is a local minimum of the path energy is 
a geodesic path, albeit not the shortest geodesic. To 
our knowledge, very few methods guarantee a globally- 
optimal solution to the problem of finding geodesics in 
shapes spaces of surfaces. Although |6i| was the first 
to provide a geometric framework for joint registration- 
comparison problem, the Riemannian metric used there 
has a limitation that it was not translation invariant. 

To handle the translation issue mentioned above, 
Jermyn et al. 121 introduced a comprehensive Rieman¬ 
nian metric that has several improvements, including the 
fact that it was translation invariant and allows some 
physical interpretations in its use. This metric, given 
later in Eqn. 0 , has terms that can be interpreted as 
measurements of bending, stretching, and changes in 
local curvatures of surfaces. (We elaborate on this topic 
later in Section |3.4[ ) It has been termed an elastic metric 
because it is invariant to reparameterizations and the 
physical interpretations associated with it. Although ||9l 
introduced this metric, it did not use the full metric 
to compute geodesic paths. Instead, it defined a new 


map, termed the square-root normal field, given by 
q{s) = ^Jr{s)nf{s) where n/(s) denotes the unit normal 
to the surface at the point 5 G S^. The square-root 
normal field has the property that the last two terms of 
the elastic metric transform to the metric under the 
map f ^ q, for some weighting of last two terms in the 
metric. The first term of the metric is discarded in this 
analysis. The transformation to metric is useful since 
one can apply some common tools from Hilbert space 
analysis to this problem, including the optimization over 
the reparameterization group for optimal registration, 
but this mapping f q is not onto and, hence, not 
invertible. The optimization step is challenging because 
the reparameterization group is an infinite-dimensional 
Frechet Lie group, and the exponential map is not a local 
diffeomorphism. Since the first term of the elastic metric 
introduced in |^l| is not used by Jermyn et al., it can result 
in zero shape distance between two surfaces that actually 
have different shapes. For example, a thin-tall cylinder 
and a fat-short cylinder, with same surface areas and 
unit normals, will have zero shape difference under this 
framework. 

Another line of work in shape analysis comes from 
Michor et al. cqi, Bauer et al. lUTl , |l2, |[T3l and Fuchs 
et al. HH (see also GSI for an overview of a lot of math¬ 
ematical results in this area). Different types of metrics 
have been studied : Sobolev metrics in CHI, curvature 
weighted metrics in tTTl , almost local metrics in [121, 
metrics mesuring the deformations of the interiors of 
shapes in (Ml . Let us mention that the first two terms of 
the metric we use in the present paper fit in the general 
study laid out in (131, and are related to the metrics 
studied in (161, (iZl, GSl (in a sense that we will make 
clear in Section |3.4l . In this set of papers, the idea is 
to replace the problem of solving the geodesic equation 
on shape space by the equivalent problem of solving 
the equation for horizontal geodesics in the pre-shape 
space. A geodesic in pre-shape space is horizontal if it is 
orthogonal to the orbits of the reparameterization group. 
One task in this strategy is therefore to compute the hor¬ 
izontal space on which the quotient map is an isometry, 
or equivalently solve a minimization problem for the in¬ 
finitesimal energy. Depending on the Riemannian metric 
on the pre-shape space, this task may be computationally 
trivial or extremely difficult to implement (for metrics 
used in im and m it is just the space of normal vector 
fields, but for metrics used in in and m it involves 
the inversion of a pseudo-differential operator). Another 
main contribution of these authors is to give sufficient 
conditions under which the Riemannian metric induced 
on shape space separates points, i.e. gives a non-zero 
geodesic distance between pairs of different shapes (a 
condition that is necessary to make shape comparison). It 
is worth noting that, in this infinite-dimensional context, 
vanishing geodesic distance is a common phenomenon 
(as was first highlighted in [TOl ). For the metric we 
use, non-vanishing geodesic distance is guaranteed by 
the non-vanishing geodesic distance on the space of 
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Riemannian metrics proved in IT^ (at least on pairs 
of shapes inducing different pull-back metrics on the 
sphere, which is what we are interested in practice). 

To summarize, the past approaches involving Rieman¬ 
nian geometry have tended to perform shape analysis 
in two steps. First, they select a representation space, 
or a pre-shape space, for objects of interest - curves 
[il, 0, GQI and surfaces El, 13, CD, ES, 113 - and 
impose a Riemannian structure on it ensuring that the 
actions of shape-preserving groups are by isometries. 
Next, they inherit this metric to the quotient space of the 
pre-shape space modulo the requisite groups, called the 
shape space, and seek geodesics between objects in this 
shape space. The task of inheriting Riemannian metrics 
to quotient spaces is complicated because reparameteri¬ 
zation groups are Frechet Lie groups and the process of 
inheriting a metric requires closed orbits, as can be seen 
in m, m, mi m, etc. Even though endowing shape 
space with a Riemannian metric (with positive distance 
function) seems to be a good approach, inducing this 
metric by a Riemannian metric on pre-shape space leads 
to difficulties that one would like to avoid (recall that 
we are only interested in shapes and not in the way they 
are parameterized). We will pursue a different strategy 
where the Riemannian metric is directly imposed on the 
quotient space, thus avoiding the need to satisfy condi¬ 
tions for inheriting metrics from the pre-shape space or 
computing an abstract horizontal space. Motivated by an 
easy implementation of the metrics, we take the point of 
view where the space of interest is the space of normal 
vector fields (in contrast with the horizontal space of a 
Riemannian submersion). Let us emphasize that there is 
no restriction in doing so : any Riemannian metric on 
shape space can by expressed as a metric defined on 
normal vector fields. 

1.2 Goals and Contributions 

Now we present the goals and contributions of this 
paper, and start by revisiting the question: What should 
be a good Riemannian metric on shape space ? A good 
Riemannian metric on shape space should be such that: 
(1) it induces a positive distance function on shape space, 
i.e. the infimum of the lengths of paths connecting two 
different shapes should be non-zero ; (2) the distance 
between two shapes should be independent of the way 
the two shapes are parameterized ; and, (3) the length 
of a path of shapes should be independent of the way 
shapes along the path are parameterized. The last point 
should be thought of as the natural generalization of the 
fact that, on a finite-dimensional Riemannian manifold, 
the length of a curve is independent of the way the curve 
is parameterized. It should be true for any path (not only 
for geodesics), and is called gauge invariance. Indeed the 
use of parameterized surfaces in order to measure the 
deformation of a shape can be compared to the use of 
a gauge. Let us comment on Fig. in order to illustrate 
this idea. Each column depicts an orbit under the repa¬ 
rameterization group for the corresponding surface, the 


surfaces in a given orbit correspond to the same shape 
but with different parameterizations. A path of shapes 
can be lifted in many ways to a path of parameterized 
surfaces. In Eig. two lifts of the bottom line path are 
depicted. The first path connects parameterized surfaces 
with different ''heights" in the fibers. This is made to 
emphasize that the variations of the "height" (i.e. of the 
parameterization) in the fibers should not influence the 
value of the length of the path of shapes. 

The main contributions of this paper are following: 

• The proposed method achieves gauge invariance, 
i.e. the lengths of paths (geodesics or otherwise) 
measured under this metric are invariant to ar¬ 
bitrary reparameterizations of shapes along these 
paths (in particular, the two paths in Eigj^have the 
same length). 

• It uses an elastic metric that accounts for any defor¬ 
mation of patches to define and compute geodesic 
paths between given objects in the shape space, and 
it presents a geometric interpretation of the different 
terms involved in this metric. 

• By defining a metric directly in the shape space, 
it avoids the optimization step over the reparam¬ 
eterization group and difficult mathematical issues 
arising from inheriting a metric from pre-shape 
space. 

Note that the third point leads to more efficient Algo¬ 
rithms in cases where one only needs a shape geodesic 
and not the optimal registration between surfaces. It 
provides the same geodesic path despite arbitrary initial 
parameterizations (or registrations) of given surfaces, 
and saves the computational cost of finding a regis¬ 
tration. This fact is also a source of limitation in the 
situation where one needs a registration. If one wants 
to use geodesic lengths for comparing shapes, then a 
registration is not needed. However, if one wants to 
study statistical summaries of deformation fields, then 
a registration will be needed. 

The rest of this paper is organized as follows. Section 2 
describes the mathematical representation of embedded 
surfaces and establishes mathematical setup. Section 3 is 
devoted to the description of gauge invariance and to 
the definition of the Riemannian metric involved in this 
paper. The geodesic computation is described in Section 
4 and Section 5 presents the experimental results. 

2 Mathematical Setup 

2.1 Notation 

We will represent a shape S with an embedding / : 

^ such that the image /(S^) is S. The function / 
is also called a parameterization of the surface S. 

We will use local coordinates {u,v) on the sphere. Eor 
the theoretical framework, any coordinates on the sphere 
are suitable, but in the application we use spherical 
coordinates : u stands for the polar angle and ranges 
from 0 to TT, and v denotes azimuthal angle and ranges 
from 0 to 27r. 
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Recall that a map / : ^ is an embedding 

when: for any point {u,v) G (1) / is smooth, in 
particular the derivatives fu and fy of / with respect 
to u and v are well-defined, ( 2 ) / is an immersion, i.e. 
the cross product fu x fy never vanishes and allows us 
to define the normal (resp. tangent) space to the surface 
/(S^) at a point f{u,v) as the subspace of which is 
generated by (resp. orthogonal to) fu^ fv, and (3) / is 
an homeomorphism onto its image, i.e. points on /(S^) 
that look close in are images of close points in S^. If 
/ is an embedding, then the surface /(S^) is naturally 
oriented by the frame {fu,fv}f or equivalently by the 
normal vector field fuXfv 

We define the space of all such surfaces as 

T '.= {f : ^ / is an embedding}. 

It is often called the pre-shape space since objects with 
same shape but different orientations or parameteriza- 
tions may correspond to different points in T. The set 
T is itself a manifold, as an open subset of the linear 
space ^^(§^, 1 ^^) of smooth functions from to 
(see Theorem 3.1 in da and the references therein). The 
tangent space to T at /, denoted by T^T, is therefore 
just C^(S2,M3). 

The shape-preserving transformations of 3D object can 
be expressed as group actions on T. The group M+ 
with multiplication operation acts on T by scaling : 
(yd, /) i-G yd/, for yd G and f ^ T. The group with 
addition as group operation acts on T, by translations : 
{v, f) f-G f V, for V G and f ^ T. The group 
SO(3) with matrix multiplication as group operation acts 
on T, by rotations : (O,/) i-> Of, for O G SO(3) and 
f ^ T. Finally, the group T := Diff^(S^) consisting of 
diffeomorphisms which preserve the orientation of 
acts also on T , by reparameterization : (7,/) / o 7“^, 

for 7 G Diff^(S^) and f e F. The use of 7“^, instead 
of 7, ensures that the action is from left and, since the 
action of SO(3) is also from left, one can form a joint 
action of G := Diff^(S^) x SO(3) x on T. In this 
paper, the translation group is taken care of by using a 
translation-independant metric (the elastic metric) and, 
when needed, the scaling is taken care of by rescaling 
the surfaces to have unit surface area. Therefore, in the 
following we will focus only on the reparameterization 
group r and on the rotation group SO(3). 

2.2 Shape Space as quotient space 

Since we are only interested in shapes of surfaces, 
we would like to identify surfaces that can be related 
through a shape-preserving transformation. This is ac¬ 
complished using the notion of group action and orbits 
under those group actions. 

Given a group G acting on T, the elements in T 
obtained by following a fix parameterized surface f 
when acted on by all elements of G is called the G- 
orbit of / or the equivalence class of / under the action 
of G, and will be denoted by [/]. In particular, when G 


is the reparameterization group, the orbit of / G is 
characterized by the surface /(S^) = S, i.e. the elements 
in [/] = 1 / 07 “^ for 7 G Tj are all possible parame- 
terizations of S. For instance in Fig. the first column 
contains some parameterized horses that are elements of 
the same orbit. The set of orbits of under a group G 
is called the quotient space and will be denoted hy TjG. 
The quotient space of interest in this paper is called shape 
space and is defined as follows. 

Definition 1: The shape space S is the set of oriented 
surfaces in which are diffeomorphic to §^, modulo 
translation and rotation. It is isomorphic to the quotient 
space of the pre-shape space by the shape-preserving 
group G := Diff+(S2) x SO(3) x • 5 = jr/Q_ 

It is important to note that the shape space S = 
J^/G is a smooth manifold and the canonical projection 
Ft : ^ J’/G, / i-G [/] is a submersion (see for 

instance 1221 and (231 ). This submersion is useful in 
establishing the notion of a vertical space that will be 
needed a little later. By definition, the vertical space 
of a submersion is the kernel space of its differential. 
When the submersion is a quotient map by a group 
action, the vertical space is the tangent space to the orbit 
(the terminology comes from the fact that the orbits are 
usually depicted as vertical fibers over a base manifold 
which is the quotient space, see Fig[^. In the case of 
the submersion II : J’/Diff^(S^), the vertical space 

takes a very natural, intuitive form. 

Proposition 2: The vertical space Ver{f) of H at some 
embedding f e is the space of vector fields which 
are tangent to the shape /(S^), or equivalently the space 
of vector fields such that the dot product with the unit 
normal vector field rif \= ^ vanishes : 

Ver{f) = {5f : ^ M^|J/(s) • nf{s) = 0,Vs G S^}. 

Remark 3: A canonical complement to this vertical 
space (consisting of tangent vector fields) is given by 
the space of vector fields normal to the surface /(S^) 
denoted by Nor. This is the sub-bundle of the tangent 
bundle TT defined by 

Nor{f) = {6f : ^ M^|(i/(s) x nf{s) = 0,Vs G S^}. 

Any tangent vector 6f G TfT admits a unique de¬ 
composition Sf = Sf'^ A Sf^ into its tangential part 
6f^ G Ver{f) and its normal part 6f^ G Nor{f). 
Specifically, the normal part is given by: 

• ( 1 ) 

See Fig. for an illustration of this decomposition. Gen¬ 
erally speaking, one has TT = Ver^Nov as a direct sum 
of smooth fiber bundles over T. This decomposition is 
preserved by the action of the reparameterization group 
r, i.e. {Sf = 6f^ o 7 and {Sf = Sf^ o 7 

(for a proof of this statement, see Section 1 of the 
Supplementary Material). 

The interest in splitting a perturbation Sf into its 
normal and vertical components comes from the fact that 
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the vertical component 6f^ G Ver{f) can only lead to 
a shape-preserving transformations of the surface /(S^). 
Thus, in the process of deforming one shape into another 
(for instance along a geodesic path) and quantifying 
shape differences between them using geodesic lengths, 
we are not interested in measuring deformations that 
are in Ver{f). An important novelty of this paper is that 
the eventual Riemannian metric is imposed only on the 
Sf^ components of the perturbations, and that the Sf^ 
components have a zero contribution to the metric. 

3 Gauge Invariance and Riemannian 
Metric 

As mentioned earlier, another important goal of this 
paper is in developing a framework that is gauge in¬ 
variant. To appreciate the utility of this framework, we 
first provide a precise definition and then motivate its 
use in shape analysis. 


3.1 Defining Gauge Invariance 

The gauge invariance relates to the parameterization 
of surfaces along a path in T and, thus, the mathe¬ 
matical objects of importance in this section are paths 
: [0,1] f-G J*. The set of such paths is the smooth 
manifold V := ^^([0, 

An element of V can be thought of as a metamorphosis 
from the initial shape to the final shape. For instance. Fig. 
[T] shows two elements in V as two different deformations 
from a parameterized horse to a parameterized cat. 
To have a picture in mind, consider the upper path 
^ : t ^{t) in V: at each time step t G [0,1], is a 
parameterized shape, i.e. a map from our model mani¬ 
fold into The map T^(0) is the parameterization of 
our initial parameterized shape chosen to be a horse and 
T^(l) is the parameterization of our final parameterized 
shape which, in this case, is a cat. 

The definition of length of the path t i-g requires 
specification of a metric on T. Given such a metric ((•,•))/ 
one can define the length as: 

L[%= ( 2 ) 

where is the velocity vector of the path 

t i-G ^(t), i.e. an infinitesimal deformation of the param¬ 
eterized shape ^(t). The geodesic distance between two 
shapes /i and /2 is then defined by 




inf im, 


(3) 


where the infimum is taken over all paths connecting 
shape /i and shape / 2 . 

We would like the length L[^], for any path to 
match the length of the path t o j(t), where 

t 7 (t) G r is any time-dependent reparameterization 
of §2 : 


= L[§], where ^(t) = o 7 (t). (4) 


More formally, set T = Diff^(§^) and define the group 
Q := G^([0, l],r), of time-dependant reparameteriza- 
tions that acts on V according to 

GxV V 

I—^ (t o7(t)). 

The group Q is called the gauge group, and one says 
that G acts by gauge transformations. We are looking for 
a framework where the length of a path is invariant to 
gauge transformations, i.e. satisfies Eqn. One should 

distinguish these transformations from temporal repa- 
rameterizations of the path ^ itself. A gauge transfor¬ 
mation changes spatial reparameterization of surfaces, 
while preserving shapes, along the path, while a tempo¬ 
ral reparameterization changes the time it takes to reach 
each shape along the path. 

To build a gauge invariant framework, the basic idea 
is as follows: take any T-invariant Riemannian metric 
((•,•)) on the pre-shape space, and ignore the direction 
tangent to the reparameterization orbit. (An example 
of T-invariant Riemannian metric is the elastic metric 
defined in Eqn. 0 as is shown in Section 2 of the 
Supplementary Material). More precisely, let ((•,•)) be a 
Riemannian metric on pre-shape space T which is pre¬ 
served by the action of the group of reparameterizations 
r, that is: 


((Vi ° 7, Sf 2 o 7 ))/o^ = {{5h,5h))f, (5) 

for any f e T, for any (i/ 1,^/2 G TfT and any 7 G T. 
Given a T-invariant sub-bundle H of TT such that 

i^(/)e^er(/)=T;^, (6) 

denote by pn ' TfJ^ H{f) the projection onto H{f) 
with respect to the direct sum decomposition given in 
Eqn. 0- This means that any element Sf G TfJ^ admits a 



Fig. 2. a. Direct sum decomposition H(f) 0 Ver(f) = TfJ^. b. Vector field 
decomposition into tangent and normal directions 


unique decomposition into the sum of an element pnidf) 
in H{f) and an element in Ver{f). We illustrate this 
decomposition of vector spaces in Fig. |^a, while the 
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particular case when H is the space of normal vector 
fields Nor is shown in Fig. |^b. 

Proposition 4: The non-negative semi-definite inner 
product on pre-shape space defined by 

satisfies the gauge-invariance condition given in Eqn. 
and induces a Riemannian metric on quotient space S 
such that the quotient map is an isometry between H{f) 
and the tangent space 


3.2 Distinction between Gauge Invariant Framework 
and Quotient Riemannian Framework 

In practice the subbundle H has to be chosen in order 
to make the implementation easy. A natural choice of 
subbundle H is the normal bundle Nor which is pre¬ 
served by the action of the reparameterization group 
r (for a proof of this statement, see Section 1 of the 
Supplementary Material). We have used this subbundle 
in the present paper. Another requirement is that the 
chosen Riemannian metric has to be T-invariant. This is 
the case for the elastic metric defined in next section. 
We will therefore apply the idea of gauge invariance 
to the concrete example of the elastic metric and the 
normal bundle Nor in the remainder of this paper. It 
is worth noting that the Riemannian metric on shape 
space obtained by restricting a Riemannian metric on 
preshape space to the normal bundle Nor differs in 
general from the quotient Riemannian metric. In fact, 
the quotient metric coincides with the restriction to the 
subbundle Nor if and only if the Horizontal subbundle 
defined by Hor{f) = Ker(d7r)-^ is the normal bundle. 
This is not the case for the elastic metric. We also remark 
that the present gauge invariant framework has been 
used implicitly in ITll , Section 6, and mil. Section 11, 
in the case where the horizontal bundle coincides with 
the normal bundle. 


3.3 Elastic Riemannian Metric 

Next, we will choose a Riemmanian metric on T that 
will enable a gauge-invariant analysis as stated above. 
We will use the elastic Riemannian metric proposed 
by Jermyn et al 151 arid given in Eqns. 0 and 
However, before we use this metric we motivate its use 
by making a connection between the space of parame¬ 
terized surfaces T and the space of metrics on a domain, 
and we will provide some geometrical interpretation of 
terms in that elastic metric. The space of positive-definite 
Riemannian metrics on will be denoted by Met(S^). 
Consider a parameterized surface / : ^ Denote 

by ^ = /*^ the pull-back of the Euclidian metric g of 
and by n/ the unit normal vector field (Gauss map) on 
5 = /(S 2 ). 


The metric g and the normal vector field n / are defined 
using derivatives of / according to: 

= (fg)^ Jac(/) = [/„ fv], and 

nf = xfj ' ||/« X fvW = Vdetff = \g\K 

where fu and fy are the derivatives of / with respect to 
the local coordinates {u,v) on the sphere. We consider 
the following relationship between parameterized sur¬ 
faces on one hand and the product space of metrics and 
normals on the other : 

^ ^ Met(S2) xC^(S2,§2) 

/ '—^ i9,nf)- 

It follows from the fundamental theorem of surface 
theory (see Bonnet's Theorem in for the local result. 
Theorem 3.8.8 in 1251 or Theorem 2.8-1 in 1261 for the 
global result) that two parameterized surfaces /i and /2 
having the same representation {g, n) differ at most by a 
translation and rotation. This is an important result, and 
implies that we can represent a surface by its induced 
metric g = f*g and the unit normal field n = rif, for 
the purpose of analyzing its shape. We will not loose 
any information about the shape of a surface / if we 
represent it by the pair {g,n). Let Sfi, 6/2 denote two 
perturbations of a surface /, and let (Sgi.Sai) = 

(^ 25 ^ 2 ) = ^*(<^/ 2 ) denote the corresponding perturba¬ 
tions in (g^n) of /. The expression for is given by: 

Sg = Jac(/)^ Jac((5/) + (Jac((5/))'^ Jac(/) 

2/« • Sfu fu ■ Sfu + fv ■ Sfu \ 

fu • Sfv + fv • Sfu 2fv -Sfv y ’ 

Sn = n +-^ ((5/„ X/„ +/„ X (5/„). 

Then, by definition, the metric on T used in the 
present paper measures these perturbations using the 
expression 

{{Sfi,Sf2})f = [ ds\g\^{aTr{g~^Sgig~^Sg2) 

+ ^ Tr( 5 “%i) Tr( 5 “^(^ 2 ) +cdni ■ (^2} . ( 7 ) 

The same metric (with a = 1 ) was introduced in 1*91, 
Eqn. (2), and called "elastic metric". A related metric 
measuring the elastic deformation of the interiors of 
shapes was used in HU (see Eqn. (4) in HU). The metric 
given in Eqn. 0 can be decomposed into three parts 

{{Sfi,Sf2}}f =f ds\g\^ {aTr{{g~'^Sgi)o{g~'^Sg2)o) 

-\-bTT{g~^Sgi) TT{g~^Sg2) + cSai • 802} , ( 8 ) 

where b = and where Aq is the traceless part of 
a 2 X 2-matrix A defined as Aq = A — ^^^^/ 2 x 2 - The 
term multiplied by a measures area-preserving changes 
in the induced metric g, the term multiplied by b mea¬ 
sures changes in the area of patches, and the last term 
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measures bending. Note that only the relative weights 
hja and c/a are meaningful. 

Now we consider a key property of this metric that 
relates to reparameterization of surfaces. Recall that 
r := Diff^(S^) denotes the subgroup of Diff(S^) con¬ 
sisting of diffeomorphisms 7 which preserve the ori¬ 
entation of i.e. such that detJac(7) > 0. (Note 
that for a diffeomorphism 7 G Diff(§^), since Jac(7) 
is invertible, the determinant of Jac(7) never vanish. 
It follows that either detJac(7)(s) > 0 for all s G 
or det Jac(7)(s) < 0 for all s G S^.) It will be called 
the group of orientation-preserving reparameterizations. 
The group T = Diff^(§^) acts on Maps(S^,]R^) by pre¬ 
composition. That is, a surface / is reparameterized by 
a 7 G Diff^(S^) according to / / o 7“^. How does 

the metric-normal representation {^g^rn) of that surface 
change due to reparameterization? This representation of 
the reparameterized surface is given by (7“^*^, n 07“^). 
This representation is T-equivariant for the actions in¬ 
troduced, i.e. if we reparameterize a surface and then 
compute its (^, n) representation, or if we compute {g,n) 
representation of a surface and then reparameterize them 
according to (7*^, n o 7), we get the same result. 

Proposition 5 : The elastic metric is invariant to the 
action of Diff^(S^). 

Proof: Please refer to Section 2 of the Supplementary 
Material. 

Although this elastic metric has been introduced by 
Jermyn et al. O, it has not been used completely for 
shape analysis of surfaces. Furthermore, we are going 
to use it in a novel way - by restricting its evaluation 
only to the normal vector fields on a surface (see next 
section for a geometric expression of the resulting metric 
on shape space). 

Definition 6: For any two perturbations Sf 1,6/2 G TfJ^ 
define the pairing 

where Sff- is the normal component of Sfi as defined in 
Eqn. Q and where ((t))/ is as given in Eqn. 

Remark 7 : It follows from proposition that ((•,•)) 
satisfies the gauge-invariant condition Lfit] = L[^], 
where is any path of shapes, ^(t) = T^(t) o with 
t 7(t) any time-dependant reparameterization, and 
Lfit] is as specified in Eqn. Q. 

3.4 Geometric expression of the eiastic metric in the 
normai direction 

In this section, we will give some geometric interpre¬ 
tation of the restriction of the elastic metric on the 
space of normal vector fields introduced in the previous 
section. Given a surface / parameterized by {u,v), we 
will consider normal variations: f^{u,v) = f{u,v) + 
eh{u,v)n{u,v), where {u,v) G §^, 5 > 0, n{u,v) = 
nf{u, v) is the unit normal to the surface /(S^) at f{u, v), 
and h : ^ M is a real function corresponding to 

the amplitude of the normal vector field hrif. Let us 


compute the first fundamental form g^ of the surface 
parameterized by /e, i.e. the metric induced on the 
parameterized surface by the Euclidian metric of 
We obtain 

fs,u '= ^ = fu^ shriu + ehuU, 

(9) 

fe,v ~ fv shriy “h shyTi. 

Therefore 


fs,u • fs,u = fu’ fu^ 2ehnu ' {h^riu • riu + hf) , 

where we have used that n • /^ = 0 and • n = 0 since 
n ' n = 1. Similarly 

f£,v ' /e,V ~ fv ' fv ^shlfly ' fy S [h Thy * Uy + , 

and 

f£,u'f£,v fu'fv~^^h {riu'fy-\-(^h Tiu'ny-\-huhy^ . 


It follows that 

g^= g^ 2sh 


'^u 

fu 


fv 

Tly 

fu 

riy 

fv 

’ 'hlu 


Tly 


' riu 

riy 

• riy 

y 




huhy A 

Khy hi ) 


Using the definition of the second fundamental form II 
of the surface /(S^), we obtain 

g^ = g- 2 ehll^e^hHlg~Hl+e^ [ hu hy )^ ( hu hy ). 


It follows that 


g-^6g = -2hg-Hl = -2hL, ( 10 ) 

where L is called the shape operator. Recall that the 
eigenvalues of L are the principal curvatures of the 
surface /(S^), denoted by and 1 ^ 2 , which provide local 
information about the surface: at a given point on the 
surface, they measure the greatest and smallest possible 
curvatures of a curve drawn on the surface passing 
through this point. For instance, the vanishing of the 
principal curvatures at one point of the surface tells that 
the surface is flat near this point (i.e. looks like a plane). 
The equality = 1^2 = l/R one point tells that the 
surface looks like a sphere of radius R near this point. 
In other words, and 1^2 are functions on the surface 
that characterize how the surface is locally curved. 

On the other hand, the variation 6 n of the normal 
vector field satisfies dn • n = 0 since the norm of n 
remains constant. Moreover n • /^ = n • /^ = 0, therefore 
6 n • fu = -n • Sfu and 5n • fy = -n • 6 fy. By Eqn. (0, 
6 fu = hriu + huTi, hence 6 n ’ fu = —hu and similarly 
6 n • fy = —hy. Consequently 6 n = afu + I6fy where 
(^) = ih^) ' follows that for two normal vector 
fields hn and kn with h,k e G^(S^,M), one has 

6 ni • 6 n 2 = K) g-^ (^- ) . ( 11 ) 

Using Eqn. jTQ| ) and Eqn. ( pT| the elastic metric restricted 
to these normal fields is given by : 

{{hn,kn))j^ = ds\g\^ \^hk (2a{Ki — K 2 )‘^ 

n 2 )‘^^c [hu hv ^ g ( 12 ) 
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This is the form used to define and compute geodesic 
paths in the shape space S in this paper. The difference 
Ki — k ,2 in the first term has been called the normal 
deformation of the surface in Ezl. The sum ki + /^2 is 
twice the mean curvature which measures variations of 
the area of local patches. These two terms are related to 
the shape index idx = -arctan ^^^^^ 12811 . The last term 
in Eqn. ( [12| ) measures variations of the normal vector 
field, i.e. bending. 

4 Geodesic Computation 

Finding geodesics between two surfaces /i and /2 under 
invariant Riemannian metrics is a difficult problem. In 
the present case, analytical solutions are not known 
and we will use a path-straightening approach to find 
geodesics. This method has been used for instance in 171 
and m. The basic idea here is to connect fi and /2 by 
any initial path and then iteratively straighten it until it 
becomes a geodesic. The update is performed using the 
gradient of an energy function. As mentioned earlier, this 
method only achieves a local minimum of the energy 
function, resulting in a geodesic path that may not be 
the shortest geodesic. 

4.1 Removing rotations and transiations 

Since we are only interested in shapes of objects and not 
in the way objects are oriented or placed in the ambient 
space we have to remove the actions of the rotation 
and translation groups. In theory, we could deal with 
them as we do with the group of reparameterizations. 
However, since SO(3) xi is just a 6 -dimensional Lie 
group in comparison to the infinite-dimensional Frechet 
Lie group Diff^(§^), it is more efficient to do the fol¬ 
lowing. First find the best translation and rotation that 
align the two objects to be compared, and then find 
the geodesic between them. To center an object we use 
Algorithm 3 given in Section 3 of the Supplementary 
Material to compute the center of mass and then subtract 
it from the surface coordinates. 

There are many ways to find the best rotation that 
aligns two shapes. In the case of elongated objects (which 
was the case in our experiments), one can do the follow¬ 
ing. Given two shapes Si and 82 ^ find the best ellipsoids 
El and E 2 that approximate the cloud of points defining 
Si and S 2 respectively, and the unitary matrices Ui and 
U 2 that map the reference axes to the axes of the el¬ 
lipsoids (with decreasing lengths). The unitary matrices 
Ui and U 2 are uniquely defined if the approximating 
ellipsoids are triaxial (i.e. the lengths of their principal 
axes are distinct). Then, we can apply the product matrix 
U 2 Ui^ on the shape Si to rotationally align with 82 - If 
one encounters a 180 degree flip, apply instead U 2 RU 2 ^, 
where R is the 180 degree rotation around the z-axis. As 
an example. Fig. shows two hands that have different 
orientations in space, the corresponding ellipsoids, and 
the hands after rotation (with a gap to separates them in 
order to facilitate visualization). 





Fig. 3. Rotational alignment: two hands before and after the alignment, 
respectively at the left and at the right. Each hand is approximated by an ellipsoid. 
The rotation used apply the axis of one ellipsoid to the axis of the other. 


To find the best ellipsoid that approximates a surface 
S and the corresponding rotation U, one can use a 
singular value decomposition of S^S. However, in the 
case where the surface is the boundary of a 3D-volume, 
it is more accurate to compute the mean of S^ S over the 
inscribed volume. It also has a more physical meaning 
since the resulting ellipsoid is equivariant with respect 
to affine transformations (see Section 3 of the Supple¬ 
mentary Material where the effect of the rotations of the 
initial surface on the ellipsoid is illustrated (Fig. 2) and 
where detailed Algorithms are presented). Moreover, the 
estimation of ellipsoid for an inscribed volume is more 
stable under reparameterizations. To illustrate this ro¬ 
bustness we show in Fig. [^different parameterizations of 
a horse (middle row) obtained by pre-composing a given 
parameterization by a diffeomorphism of the sphere 
(bottom row) and the resulting ellipsoid (top row). The 
diffeormorphims used in this experiment are (from left 
to right) (pi = identity, ip 2 = rotation of — 37 r /4 around 
x-axis, (p 3 = Mobius transformation that maps 2 : G - 
C U {oc} to = 0.4^ + 0.5, (p 4 = rotation of — 7 r /2 

around x-axis composed with (^ 3 . 



Fig. 4. Robustness of the approximating ellipsoid of a surface with respect to 
reparameterizations. 


In the case where the approximating ellipsoids are 
not triaxial, one has to use additional information about 
the surfaces to align them properly (for instance, one 
can use four points on each surface). This case was not 
implemented in the present paper. 

4.2 Computations of the energy 

Let : [ 0 , 1 ] ^ R. The energy of the path is defined 
to be: 

dt = ((^i, dt, 
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where ((•,•)) is the elastic metric given in Eqn. 

’ 'ri)n is the normal component of the de¬ 
formation, and ((•,')) is the inner product presented in 
Eqn. We will present several numerical strategies 
for approximating this energy and will compare their 
computational costs in Table This evaluation uses a 
linear path connecting two concentric spheres of radius 
Ri = 1 and R 2 = 2.5, with constants a = 1, A = 0.125 
and c = 0 for defining energy (see Eig. |^. The the¬ 
oretical value of the energy in this case is given by 
Eth = 327r(a + A)(i ?2 — Ri)‘^ and measures exclusively 
the cost of changing the area of the spheres (the first and 
third term of the metric given in Eqn. ^ vanish in this 
experiment). We expect that improvement in accuracy 
comes at an increased computational cost, and this is 
indeed the case in the results presented in the Table. Note 
that a time-dependent rotation is applied on the path of 
spheres, but the values of the energy is independant of 
this rotation. 



Fig. 5. Path connecting two concentric spheres used for computations in Table 

ffl 

One way to compute the energy of a path of shapes 
is to express it using the coefficients of the first funda¬ 
mental form. Consider the mapping ^ : §2 X R ^ m 3 
and define 

E = F = G = (13) 

and their time derivatives 

E = 2VI/4 ■^u,F = Hu H, G = 2H • 

as well as the unit normal field n := Uf = and 

the vector field w = + ^u x Then, the 

energy of a path ^ decomposes into the sum of four 
terms: f (T^(t)) = Ei E S 2 E S 3 E S 4 , where 

Si= afX^iEG-FH^GBdududt 

with B = G^E^ + 2{EG+F^)F^+E^G^ 
-AFGEF+2F^EG-4:EFFG , 

£ 2 = {Hi)lX^iEG-FHHGE-2FF+EG)^ckdvdt , 
^ 3 = -c^l y{GE-2FF+ EG){n-w)dadvdt , 
cfl H(EG-FH^i.ww)dudvdt . 

In the implementation of these formulas, we can reach 
singularities on the boundary of the integration domain, 
which we can ignore. In the example involving two 
concentric spheres, the total energy computed by this 
method is labelled Ejg^jj in Table 

Another way to compute the energy is based on 
Eqn. that expresses the elastic metric in terms of 
principal curvatures. In terms of the coefficients of the 
first fundamental form given in Eqn. ( [T^ and of the 


second fundamental given by 

^ ~ ^UU ’ 

f = ^uv • n = -^u • riy = • Uu, 

g = ^yy • n = -^y * Uy , 

the Gauss curvature K and the mean curvature H have 
the following expressions 

eg-f l eG + gE-2fF 

EG-F^’ 2 EG-F^ ’ 

and the principal curvatures are given by 

ki=H+ - K, K2 = H- - K. 

Again, in the implementation of these formulas, we can 
get singularities for curvatures on the boundary, but 
we can ignore them in computing the integral given 
in Eqn. This corresponds to removing a small disc 
on the parameterized surface around the images of the 
north and south poles. 

In the example of the two concentric spheres, the 
theoretical values of and H 2 is the constant function 
equal io 1/R where R = Ri H t{R 2 — Ri) is the radius 
of the sphere along the path interpolating linearly the 
sphere of radius i?i = 1 to the sphere of radius i ?2 = 2.5. 
The total energy computed by this method is labelled 
in Table ^ 


Fig. 6. From left to right: A hand with the tangent plane and normal at the tip 
of the index finger; 3-neighborhood of the tip of the index finger; tip of the index 
finger after rotation; a closeup; approximating second order polynomial. 

To improve the computation of the curvatures and 
therefore also of the energy, we can use polynomial ap¬ 
proximations of the surfaces. This procedure, leading to 
the computation of the principal curvatures, is illustrated 
in Eig. To compute the principal curvatures at a given 
point of a surface, e.g. at the tip of the index finger 
of the hand depicted in Eig. we first compute the 
normal at this point by averaging the normals of the 
facets having this point as vertex. A tangent plane is then 
defined as the plane orthogonal to the normal passing 
through the point under consideration. A neighborhood 
of the point is isolated from the surface (we use a 3- 
neighborhood, see second drawing in Eig. |^. We then 
apply a rigid transformation to center the point at the 
origin and to align the tangent plan with the xy-plane 
(see third drawing, and a closeup in the fourth drawing). 
After that, we use Algorithm 5 given in Section 5 of the 
Supplementary Material to compute the second order 
polynomial P(x, y) = aix^ + a 2 y^ + ctsxy + a 4 X + a^y + ae, 
which minimizes the sum — P{xi^yi)Y over the 

points of the centered and rotated neighborhood. Then, 
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the Gauss curvature at that point is given hy K = 
4 aia 2 — ag, the mean curvature by = ai + a 2 , and the 
principal curvatures by ni = ai + a 2 + ~ <^ 2 )^ + <^ 3 ) 

and />:2 = + a 2 — — <^ 2 )^ + o^)- In the example of 

the two concentric spheres, the total energy computed 
using the principal curvatures obtained by this method 
is labelled Ep in Table [T] 



Fig. 7. A path of zero energy connecting a hand and the same hand with 
another parameterization. 

In order to show that the energy function of a path 
of shapes is independent of the way the objects are 
parameterized, we replace the integration over the do¬ 
main of parameterization by the integration over the 
triangulated surfaces. This means that we approximate 
the area elements of the surfaces by the area of triangles 
whose vertices are given by the parameterization. In 
this way, the parameterization of surfaces is only used 
to define the surfaces, but plays no role at all in the 
computation of the energy function. In the example of 
the two concentric spheres, the total energy computed 
by this method is labelled E/\ and given in Table 
In Fig. 0 a path connecting a hand to the same hand, 
but with a different parameterization, is shown. The 
energy of this path, computed with the constants a = 1, 
A = c = 0.125, reads = 0.4824, hence is close to 
0. Now returning to Fig. the energy of the lower 
path from a horse to a cat computed with the same 
constants is = 227.4049, its length is L[^] = 14.9099, 
whereas the upper path (obtained from the lower path 
by applying a different reparameterization at each time 
step) has an energy equal to E/\ = 225.5249 and a length 
of L[^] = 14.8802. Note that in this example, the colors 
refer to the Euclidean distance to the point on the surface 
corresponding to the image of the north pole of the 
sphere (cold colors for small distances versus hot colors 
for large distances). In particular, the north and south 
poles do not correspond in these two paths. 

4.3 Orthonormal Basis of Deformations 

In this section, we define bases for representing per¬ 
turbations of a path of surfaces. These basis elements 
form possible directions for use in path-straightening in 
Section 14.41 The first basis we used is a variation of the 
one given in m. We start with a basis Si = {Yj^, 1 < 
I < N,—I < m < 1} oi spherical harmonics of degree less 
than N, available in Matlab as function SPHARM (see 


Energy, 10"^ points per object 

Elapsed time for 10^ points 

Emu = 246.2854 

0.221726 seconds 

= 249.1969 

0.862376 seconds 

Ep = 255.8288 

1.238354 seconds 

Ea = 255.9043 

9.738431 seconds 

Energy for 4 x 10"^ points 

Elapsed time for 4 x 10"^ points 

Emu = 249.1503 

0.978828 seconds 

Ek,k2 = 251.8494 

3.45599 seconds 

Ep = 254.7646 

4.906798 seconds 

Ea = 254.7832 

39.011899 seconds 


TABLE 1 

Computation of the energy of a path connecting two concentric spheres (Fig.|^ 
using different methods, and time needed for the computations. The theoretical 
value of the energy is Eth = 254.4690. Here Ri = 1, R 2 = 2.5, A = 0.125 
and c = 0. 


|[29l for more information on spherical harmonics). We 
make three copies of this basis of M-valued functions 
in order to obtain a basis S 2 of the space I/^(§^,M^) 
of -valued functions. Similar to Xie et al. ISOl , we 
demonstrate reconstruction of some surfaces using the 
resulting basis, as the degree of the spherical harmonics 
grows, in the Supplementary Material (Fig. 5). 

Next, we want to construct a basis of perturbations 
of a path connecting two parameterized surfaces fi and 
/ 2 . In order to apply the path-straightening method as 
described in Section |4.4} we want the perturbations to 
vanish at t = 0 and t = 1 so that fi and /2 remain 
fixed. Therefore, we want a basis of I/^(S^ x [ 0 , 1 ],M^) 
with elements that have this property. To ensure this, 
each element of 82 is multiplied by a basis element of 
I/^([ 0 , 1 ],M) of the form Pj{t) = ^ sin( 7 rjft), I < j < J. 
Unfortunately a major limitation of the resulting basis 
is that slowly- and rapidly-oscillating harmonics have 
comparable amplitudes. In the implementation of the 
path-straightening method, this implies that the updated 
path can go out of the open set of immersions. 

One possible way to counter this effect is to or- 
thonormalize the -basis with respect to an i7^-type 
scalar product (i.e. that measures also the variation of 
the derivatives). For this kind of scalar product, an 
orthonormal basis consists of functions which have con¬ 
trolled derivatives (hence can not oscillate to much). 
This approach was also used in |7l where the L^-basis 
is orthonormalized using the Gram-Schmidt procedure 
with respect to the following scalar product 

{B\B^) = /; L (B^ ■ B^ + Bl ■ B^ + BI-BI 

+ ^t,u ' ^i,u + Bt,v ’ ^t,v) 

However, when increasing the degree of spherical har¬ 
monics, the computational cost of generation of an or¬ 
thonormal basis using this scalar product is very high. 
Therefore, we first orthonormalize the basis 82 with 
respect to the following inner product 

{B\B^) = /g, (B^-B^+Bl-Bl + Bl-B^,) ds, (14) 
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and then we multiply the resulting basis by the time- 
dependant components Pj{t) = | sin( 7 rjt), I < j < J. 
The advantage of this method is that the Gram-Schmidt 
procedure is applied to matrices of lower dimensions 
(without the time dimension) and on a smaller number 
of elements (by a factor J). The spatial oscillations of 
the resulting basis elements are well controlled by the 
presence of the spatial derivatives and By in the inner 
product given in Eqn. 


4.4 Path-straightening method 

The path-straightening method is used to find critical 
points of the energy functional. Starting with an arbitrary 
path, the method consists of iteratively deforming (or 
"straightening") the path in the opposite direction of 
the gradient, until the path converges to a geodesic. 
The gradient of the path energy is approximated using 
a basis B of possible perturbations of a path of sur¬ 
faces as constructed in the previous section. We first 
compute the directional derivatives = J^(f (T^ + 

eb))\^=o where b ranges over B. This is done by fixing 
a small ei and approximating the directional derivative 
by Vf^( 6 ) ^ {£{^ + eib) - £{^))ei~^. Using the finite 
orthonormal basis S, we obtain a numerical approxima¬ 
tion of the gradient: Vf^( 6 ) b. In particular, 

the norm of the gradient is approximately given by 
llVf^lP = update of the path is 

done by replacing ^ by — e 2 V 5 ^, where 62 is a small 
parameter that has to be ajusted empirically. The method 
is detailed in Algorithm below. 


Input: 

1) A path ^ between two parameterized surfaces /i and / 2 , 

2) a basis of perturbation B. 

Output: 

1) The minimal energy needed to deform fi into /2 given by the value of 
the cost function E, 

2) the geodesic path between fi and / 2 . 

Set \\VEf = 1. 

while ||V£;||^ > 10“^ do 

2- Compute the energy E of the path ^ according to Eqn.j^ or 
Eqn. 

3- Set ^upd = 0 and || = 0. 

for i ^ 1 to size{B) do 

4- Add a perturbation to the current path define 

qf(i) = ^(*)/ where B(i) is the element of the 

perturbation basis B of index i and ei > 0 is small. 

5- Compute the energy E(i) of the perturbed path ^(i). 

6- Compute the gradient of energy VE(i) in the direction B(i) 
using the approximation VE(i) ~ e{z)-e ^ 

7- Compute the updating path: ^ ^upd E VE(i) ■ B(i). 

8- Compute the squarred norm of the gradient of energy at path 
'I': IIVBf ^ IIVEf + (VE{i)y. 

end 

10 - Update the path: ^ ^ — e2^upd 

end 

Algorithm 1: Path-straightening method. 


5 Experimental results 

The 3D realistic models used in our experiments are part 
of the TOSCA 1311 dataset. Their spherical parameteri- 
zations were initially implemented in l32l . 


5.1 Examples of geodesics obtained by path¬ 
straightening 

First we apply the path-straightening method to the case 
where the surfaces at the extremes of the initial path have 
the same shape, but different parameterizations. More 
precisely, we consider the special case where ^o(O) = /i, 
^o(l) = fi o 7 for some diffeomorphism 7 and where 
we initialize the path with piecewise linear interpolation 
to a different surface /s in the middle of the path, i.e. 
^o(|) = /s- This situation is illustrated in Fig. The 
proposed gauge-invariant approach is expected to reach 
a path with constant shape as a geodesic, despite the 
different shapes appearing in the initial path and the 
different parameterization of shapes at the end points of 
the path (to emphasize the differences in parameteriza¬ 
tion, zoom-ins of these surfaces are also shown). Once 
we have the geodesic path ^ between the given surfaces, 
the distance in the shape space between fi and fi o 7, 
12 )f is the length of as specified in Fqn. As 
expected, the resulting geodesic path, shown in Fig. 
is constant with the same shape as the either end, and 
with f 2 ) = 0. Using path-straightening, we obtain 

a 99.28% decrease in the energy function from the initial 
path to the final path. 

In Fig. we consider more challenging shapes. The 
top-two rows display the case where we have ^o(O) = 
/i, T^o(l) = fi (a cat) and where we initialize the path 
with piecewise linear interpolation to a horse in the 
middle of the path. The upper row shows the initial 
path and the second row the geodesic path. We can see 
that the geodesic path has a constant shape throughout, 
as expected. We also plot the evolution of the path 
energy on the right during path-straightening. We can 
see that the energy decreases until it reaches a relatively 
small value; the theoretical minimum is, of course, zero 
for a contant path. In the last two rows of Fig. we 
consider the case of two hands. We initialize the path 
with linear interpolation (third row in Fig. 0, and the 
resulting path is shown in the last rows of Fig. The 
energy evolution is shown on the right and we can see 
the energy decreasing until it reaches a constant value; 
thus, the final path is a geodesic. It can be seen that 
the deformation along the geodesic path is more natural 
than the original path. 



»»»A 


Fig. 8. Illustration of initial path (upper row) and geodesic path in shape space 
(middle row). The energy is reported in the button row. The surfaces at the end 
points of the path have different parameterizations. 
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Fig. 9. The top row shows an initial path formed by linear interpolation between a cat to a horse and back to the cat. The second row illustrates the geodesic obtained 
after 800 iterations of path-straightening. The corresponding evolution of the energy is shown on the right. Similarly, the third row shows a linear path between two hands 
with bad correspondence and the last row shows the final geodesic, with the corresponding energy is shown on the right. 


5.2 Classification of 3D shapes 

As mentioned earlier, the geodesic paths provide us 
with tools for comparing, and deforming parameterized 
surfaces. We suggest a comparison of shapes of 3D 
objects using geodesic distances between their boundary 
surfaces in the shape space. This section presents a spe¬ 
cific application to illustrate that idea. In this section, we 
study several shapes belonging to four classes: horses, 
hands, cats and centaurs. 

We begin by computing the pairwise geodesic dis¬ 
tances between corresponding 3D surfaces. The distance 
matrix and the classification dendrogram are shown in 
Fig. In the distance matrix, we can easily distinguish 
four classes corresponding to four blue boxes. Actually 
the cold colors in the illustrated matrix correspond to 
small values of distances versus hot colors that corre¬ 
spond to greater distances. The clustering obtained using 
the dendrogram (command in matlab) can be interpreted 
by slicing the top of the dendrogram by a horizontal line 
to split the shapes into the desired number of classes, 
and then sliding the horizontal line to the bottom in or¬ 
der to refine the classification. The coarsest classification 
results by slicing the dendrogram into two classes (by a 
horizontal line close to the top), the shapes 4, 5 and 6 
(the hands) forms a first class and the remaining (horses, 
cats and centaurs) are grouped together as a second 
class. The next level in classification distinguishes the 
shapes 1, 2, and 3 (the horses) and 12, 13 (the centaurs) 
from the shapes 7, 8, 9, 10, 11 (the cats). The finest 
level separates the horses and the centaurs in different 
classes and results in four classes. Thus, we argue that 
the proposed framework provides a powerful tool for 
shape classification. 



Fig. 10. Classification performance; left: the distance matrix, right: the 
dendrogram. 


5.3 The effect of number of basis elements 

In this section, we study the effect of the number of basis 
elements, used in path-straightening, on the resulting 
geodesic path. Given two parameterized surfaces fi and 
/ 2 , we again initialize the path with the linear inter¬ 
polation to a different surface /s in the middle of the 
path. This initial path is shown in the upper row of Fig 
pdj Then, we compute the geodesic path using different 
number of basis elements. We show the geodesic paths 
that use 52, 432 and 1728 basis elements, respectively. 
We can see that the larger the number of basis elements, 
the better the final result is. We also provide the trade-off 
between the number of basis elements and the minimum 
energy value obtained. The trade-off confirms our asser¬ 
tion. At the bottom of the figure, we show the geodesic 
path obtained when the path-straightening Algorithm 
is initialized with the linear interpolation between /i 
and / 2 . This path is also calculated using the number of 
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basis elements corresponding to the lowest energy. This 
path can be seen as ground truth to visually interpret 
the previous geodesics (with more complicated initial 
conditions and fewer basis elements). 


pt ^ 0 ^ ^ 



200 400 600 800 1000 1200 1400 1600 1800 

# basis elements 


Fig. 11. The effect of the number of basis elements, (1) initial path, (2) geodesic 
path using 52 basis elements, (3) geodesic path using 432 basis elements, (4) 
geodesic path using 1728 basis elements, (5) geodesic path using 1728 basis 
elements after linear interpolation initialization. 


6 Conclusion 

In this paper we have proposed a novel Rieman- 
nian framework for computing geodesic paths between 
shapes of parameterized surfaces. These geodesics are 
invariant to rigid motion, scaling and most importantly 
reparameterization of individual surfaces. The novelty 
lies in defining a Riemannian metric directly on the 
quotient (shape) space, rather than inheriting it from pre¬ 
shape space, and in using it to formulate a path energy 
that measures only the normal components of velocities 
along the path. The geodesic computation is based on 
a path-straightening technique that iteratively corrects 
paths between surfaces until geodesics are achieved. We 
have presented some examples of geodesics between 
surfaces in shape spaces and utilized the distances be¬ 
tween surfaces for classification of some 3D shapes. 
However, the computational costs of our programs are 
deemed high and convergence should be accelerated in 
order to be able to apply this framework in realistic 
practical scenarios such as, for instance, human body 
action recognition. 
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SUPPLEMENTARY MATERIAL 


The Subbundle of normal vector fields 

IS r-INVARIANT 


In this section we provide the proof of the fact that 
the subbundle of normal vector fields is a F-invariant 
complement to the subbundle of tangent vector fields. 

Proposition 8: Denote by Nor the subbundle of the 
tangent bundle TT consisting of normal vector fields 
which is the space of vector fields such that the cross 
product with the normal rif : ^ to the shape 

/(S^) vanishes : 

Nor{f) = {6f : Ef ^ such that 6f x rif = 0}. 


Any tangent vector Sf e T^T admits a unique decom¬ 
position 

Sf = sf + sf 

into its tangential part 6f^ G Ver{f) and its normal part 
6f^ G Nor{f). In other words one as 

TT = Ver 0 Nor 


as a direct sum of smooth fiber bundles over J^. More¬ 
over this decomposition is preserved by the action of the 
re-parametrization group F, i.e. o 7 )^ = 6f^ o 7 and 
(< 5 / 07 )-^ = 5 /-^ 07 . 

Proof: The uniqueness of the decomposition into 
tangential and normal direction comes from the unique¬ 
ness of the decomposition of a vector in into a 
tangent vector and normal vector to the surface. The 
smoothness of the decomposition is a consequence of 
the smoothness of the tangent and normal bundles. To 
see that F preserves the normal bundle, note that if 7 G F 
reads 7 = (71(1^, 1;), 72(1^,'i;)) in a chart, then (/ o 7)^ = 
^ and (/o 7 )„ = /„07 |a +/„07 

therefore 


(/°7)„x(/°7)„ 




/ 971 ^72 
7 du dv 


dl2 dll ^ 
du dv ) 


It follows that the unit normal vector field to the 
parametrized surface / o 7 reads 

U °l)u X (/° 7 )„ 

II X II 


f ^71 ^72 ^72 dll \ 

/u o 7 X /V o 7 \du dv du dv J 

“II/.O7X/.07II' 

where in the last equality we have used that 7 preserves 
the orientation of S^. Therefore Sfoj = (^/^ 0 ^/^)o 7 = 
Sf^oj^df^ 07 with 5f^ oj G Ver{f oj) and 5f^ 07 g 
Nor{f o 7). The uniqueness of the decomposition then 
implies 6f^ 07= [Sf o 7)^ and {Sf o 7)“^ = Sf^ o 7. 


Proof of f-invariance of the elastic met¬ 
ric 

Now we will prove the fact that the elastic metric 
is invariant by the group of orientation-preserving re- 
parametrizations F = Diff^(S^). This means that 

= {{h,k))f. 

for 7 G r and any tangent vectors h^k at f e T. 

Denote by / := /07. Set {g, Uf) := $(/), and ( 5 , 7 ) := 
<F(/). Define h := h o j and k = k o j. Let us compute 
the volume form of the metric g. For any 5 G S^, one has 

Jac/(5) = Jac(/o7)(5) = {Jac f){'y{s)) • Jac7(5), and 
g{s) = (Jac7)^^(7(5)) (Jac 7 ). Therefore 

det g{s) = det(Jac 7 )^ det g{j{s)) det Jacy 
= (det Jac 7 )^ det ^( 7 ( 5 )), 


and 


| 5 (s )|2 = 7 det 5 (s) = 7 (det Jac 7)2 det g{^{s)) 

= |det Jac 7 | |5(7(s))|5- 

Let us now compute the first two terms of the elastic 
metric. Since 

g{sf = (Jac 7 )“^ g{l{s)f (Jacf) ^, 


and 


one has 


(^(s) = (Jac 7 )^(^( 7 (s)) Jac 7 , 


Trfl = 


Tr 


Jac7 ^ 5 ( 7 ( 5 )) ^ (Jac7 (Jac7) (%r(7(s)) Jac7 
= Tr (Jac7)“7(7(s))”^'%'(7('S)) Jac7 
= Tr5f-i<^(7(s)). 


Therefore, if one denotes by (resp. (^ 2 ,^ 2 )) 

the infinitesimal variation of the pull-back metric g and 
the normal vector field n /07 induced by the tangent 
vector h G Tfo^T (resp. k G Tfo^T), and (^1,^1) 
(resp. (^ 2 ,^ 2 )) the infinitesimal variation of the pull¬ 
back metric g and the normal vector field nf induced by 
the tangent vector h eTfJ^ (resp. k G T/ 7 ^), one has 


Trg ^<^ 2 ( 5 ) = 


Tr 


7ac7) 7 7i (Jac7) (Jac7) 7 72Jac7 
Tr5“7ifi'“^<%'2(7(s))- 


For the last term of the metric, since 7 acts by re- 
parametrization on the normal vector field, one has 
^ 1 ( 5 ) = ^ 1 ( 7 ( 5 )) and Sri 2 {s) = Sai{j{s)). The invariance 
by re-parametrization of the elastic metrics then follows 
by a simple change of variables in the integral defining 
it. 
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Alignment of 3D-shapes 


In this section we provide the details of the alignement 
program (Algorithm described in section 4.1 of the 
paper. The computation of the inscribed volume in a 
surface / is made using Algorithm The center of mass 
of the inscribed volume in a surface / is computed using 
Algorithm]^ The computation of the second moments is 
implemented using Algorithm More pictures illustrat¬ 
ing the robustness of the approximating ellipsoid when 
the parametrization of the initial surface is changed are 
given in Fig. The dependance of the approximating 
ellipsoid with respect to a rotation of the surface is 
illustrated in Figjl^ 





% 



^ ^ ^ 

% 


^ \ % 
\ \ 







Fig. 12. Lower rows: different re-parametrizations of the 
sphere; Middle rows: corresponding re-parametrizations 
of a shape; Upper rows: corresponding approximating 
ellipsoids. 


Input: 

1) a grid of n x n points on the unit sphere, i.e. for each index 
(z, j) G [1, n] X [1, n], a value of polar angle 0(i^j) and of 
azimuthal angle 

2) a parametrized surface /i, i.e. for each index 

(z,j) G [1,^] X [l,n], a point /i(^, j) in corresponding to 
the image of the point on the sphere with spherical coordinates 
(e(ij),(j)(ij)) by the map /i, 

3) a parametrized surface / 2 , i.e. for each index 

{i,j) G [l,n] X [l,n], a point / 2 (f, j) in corresponding to 
the image of the point on the sphere with spherical coordinates 
by the map / 2 . 

Output: 

1) a centered and scaled surface Fi having the same shape as /i, 
with center of mass at the origin and inscribed volume 1, 

2) a centered, scaled and rotated surface F 2 having the same 
shape as / 2 , with center of mass at the origin, inscribed volume 
1 and principal axes aligned with the principal axes of Fi . 

3) For A: = 1, 2, an approximating ellipse Ej. of Fj.. 

Algorithm: 

1- For A: = 1, 2, use algorithm to compute the volume Volk 
inscribed in the surface fk- 

2 - Forfc = l,2, /fc^/fc/(Voyi/3 

3- For A: = 1, 2, use algorithm 0 to compute the center of mass 
Center^ of the inscribed volume in surface fk- 

4- For k = 1,2, fk ^ fk - Center^. 

5- For A: = 1, 2, use algorithm to compute the second moments 
Mk of surface fk - 

6- For A: = 1, 2, compute [Uk, Sk, Vk] = svd(Mk). 

7- Set Fi = fi and F 2 = U 2 x U[ x / 2 . 

8- For A: = 1, 2, compute 

/ TtT \ 5 1 , - 

— 1 det(Mfc) loUkX^/^xU'k- 

9- Set Ek = Ak X sphere. A: = 1, 2. 

Algorithm 2: Alignement of 3D-shapes 


Input: 3D-parametrized surface / of size a x b x 3. 

Output: 

1) Inscribed volume Vol in surface / 

2) volume voll(i, j) of infinitesimal tetrahedron with vertices 

0 , f{i + l,i), f{i,j + 1 ); 

3) volume vol2(z, j) of infinitesimal tetrahedron with vertices 

0 , f{i -I- 1 , i -I- 1 ,:), f{i + l,i), fii, i -I- 1 ). 

Algorithm: Initialize Vol = 0. 

for z ^ 1 to size{f, 1) do 

for j 1 to size{f, 2) do 

1- Set 

edge(l) = + 

edge(2) = f(i,j -I- 1,:) - f :) 
edge(3) = f {i,j -|- 1,:) - f{i + l,j + 1,:) 
edge(4) = /(i + 1, j, :) - f(i + lj -|- 1,:) 

2- Set 

voll(i,i) = |Det(edge(l),edge(2),-/(i,i,:)) 

vol2(i,j) = |Det(edge(3),edge(4),-/(i+l,i + l,:)) 

3- Vol ■(- Vol + vol(l) + vol(2). 

end 

end 


Algorithm 3: Computation of inscribed volume 
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Input: 

1) 3D-parametrized surface / of size a x b x 3 

2) Inscribed volume Vol in surface / 

3) volume voll(i, j) of infinitesimal tetrahedron with vertices 

0, /(*, i, 0> /(* +i). f(i, j + 1) 

4) volume vol2(z, j) of infinitesimal tetrahedron with vertices 

0, f(i + l,i + 1, 0, /(* + 1, i), f(i,j + !)■ 

Output: Center of mass of inscribed volume in surface /. 

Algorithm: Initialize Center = (0, 0, 0). for z ^ 1 to size{f, 1) do 
for j ^ 1 to size{f, 2) do 

1- ml = + f{i + l,j,-) + f{i,j + 1 ,-)) 

2- m2= |(/(i + l,i +1,:) + /(! +l,i,:) + /(i,i + l,:)) 

3- Center ^ Center + voll(f, j) xml+vol2(i, j) xm2 

end 

end 

Algorithm 4: Computation of center of mass 



Input: 

1) 3D-parametrized surface / of size a x b x 3 

2) Inscribed volume Vol in surface / 

3) volume voll(i, j) of infinitesimal tetrahedron with vertices 

+ fii, j + 1) 

4) volume vol2{i,j) of infinitesimal tetrahedron with vertices 

0, f{i + l,i + 1, :)> /(* + l>i)> /(*>i + !)■ 

Output: second moments of surface / defined as the following 
integral over the inscribed volume 

/ / xy xz \ 

I xy yz dvol 

y xz xy z‘^ J 


Fig. 13. Dependance of the approximating ellipsoid with 
respect to rotation of the shape. 

Second order approximation of a sur¬ 
face 

Algorithm gives the second order approximation of a 
surface at a given point. It was used in section 4.2 in 
order to compute the principal curvatures of a surface. 


Algorithm: Initialize M = zeros(3, 3). 
for z ^ 1 to size{f, 1) do 

for j 1 to size{f, 2) do 

for ^ 1 to 3 do 
for / ^ 1 to 3 do 

= (f(i, j, k) + f(i, j + l, j, l) + f(i, j + l, l)), 


s2 = k) + f(i + l,j, k))*ifii,j, 1) + f(i + l,j, 1)), 


S3 = (/(i+1, J, k) + fii,j + l, k))*ifii+l,j, l) + f(i,j + l, 0 ). 
1 

ml = — * (si + s2 + s3). 
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s4 = (/(i+1, J, k) + fii,j + l, k))*{fii+l,j, l) + f(i,j + l, 0). 


s5 = (/(-i+l, j+1, ky\-f(i+l, j, fc)>(/(-H-l, j+1, l) + f(i+l, j, 1)), 


s6 = i+l, ky\-f(i, j+1, fc)H/(i+l, j+1, l) + f(i, j+1, 1)) , 

1 

m2 = — =)= (s4 + s5 + s6). 
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M(k, 1) ^ M{k, 1) + voll(i, j). ml + vol2(i, j). m2. 


end 


end 


end 


end 


Algorithm 5: Computation of second moments 


Input: a surface passing through the origin, tangent to the 
fcy-plane at the origin 

Output: coefficients A = (ai, a 2 , U 3 , U4, U 5 , ue) of the second 
order polynomial 

P{x, y) = aix‘^ + a 2 y‘^ + a^xy + a^x + a^y + ae, which 
minimize the sum — P{xi^yi))^ over the points 

{xi^yi, Zi) of the surface 

Algorithm: Initialize zb = zeros(l, 6), B = zeros(l, 6), 

B2 = zeros(6, 6). 

for z ^ 1 to number of points do 

1 - zb{ 1) <r- zb{1) + Zix‘j^-, 

^s(2) ^ zb{‘2) + ziy ^; 
zb{‘^) ^ zb{3) + ZiXiyi', 

Zb{^) ^ Zb{^) + ZiXi] 
zb{^) ^ zb{3) + ziyi] 

Zb{^) ^ Zb{^) + Zi’^ 

2- B{l)^xl- 
B(2)^yh 

B{3) +- Xiyi] 

B(4) ^ Xi) 

^( 5 ) Vi: 

B(6) ^ 1; 

3- B2^ B2 + B' ^ B] 

end 


A = irw{B 2 ) * z'^. 

Algorithm 6: Computation of second order approxima¬ 
tion of a surface at a given point. 
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INDEPENDANCE OF THE ENERGY FUNCTION 
WITH RESPECT TO REPARAMETRIZATION 

From the theory it is clear that the energy function 
defined in the paper is independant of the way shapes 
are parametrized along a path. To provide numerical ex¬ 
amples to illustrate this fact was however a difficult task. 
One has to mention here that the parametrization has to 
be changed smoothly in order to provide a smooth path 
in the pre-shape space. Therefore only parametrizations 
that are closed to the initial parametrization can be used 
in these experiments. 

In Figure we give examples of zero-energy paths, 
projecting to a point in shape space. The energy, as 
computed by our program, is closed to zero for each 
of them. 

In Figure we are interested in two different lifts of 
the same path in Shape space. The rows go by pairs: the 
upper two rows show a metamorphosis from a horse to 
a jumping cat, but with two different parametrizations. 
Theoretically the energy of the two upper paths should 
be the same. Numerically we obtain an energy = 
225.3565 for the upper path, and Ea = 225.3216 for the 
second one. For the third and forth paths, showing a 
metamorphosis from a jumping cat to a standing cat, 
the energy computed by our program is Ea = 180.8444 
and Ea = 176.8673 respectively. For the fifth and sixth 
paths, from a standing cat to a standing horse, the com¬ 
puted energies are Ea = 243.1812 and Ea = 239.5410 
respectively. These energies were computed with the 
parameters a = 1, A = 0.125, c = 0, 50^ numbers of points 
and using 6-neighboordhoods for the computation of 
principal curvatures. 
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Fig. 14. Four Paths connecting the same shape but with a parametrization depending smoothly on time. The energy 
computed by our program is respectively £'a = 0 for the path of hands, Ea = 0.1113 for the path of horses, £'a = 0 
for the path of cats, and Ea = 0.0014 for the path of Centaurs. 
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Fig. 15. Pairs of paths projecting to the same path in Shape space, but with different parametrizations. The energies 
of these paths, as computed by our program, are respectively (from the upper row to the lower row): = 225.3565, 

= 225.3216, E^ = 180.8444, E^^ = 176.8673, E^ = 243.1812 and £’a = 239.5410. 
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Fig. 16. Reconstruction of several surfaces with different degree of spherical harmonics. From left to right is depicted 
the initial surface and its approximation using spherical harmonics with maximal degree f = 3, f = 5, f = 7, f = 11, 
f = 15, f = 18, f = 20 and I = 28 respectively. 
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TABLE 2 

List of symbols and their definitions used in this paper. 


Symbol 

Definition /Explanation 


the unit sphere in 

(u,v) 

local coordinates on 

a • b 

the Euclidean inner dot product in 

a X b 

the cross product in 

S 

an observed surface modulo translation and rotation, or shape. 

s 

the set of observed surfaces in modulo translation and rotation or shape space. 


length of the path of surfaces. 

distlSi,S 2 ) 

infimum of the lengths of paths in S between S\ and S 2 , 

TsS 

the tangent space to the space of shapes <S at a given shape S, 

It is identified with the space of all vector fields orthogonal to the surface S, 

{XuX2) 

a Riemannian metric on the shape space S, 

MSi,S2 ) 

a geodesic path in <S, from to 82 , parameterized by t G [0,1], 

f or F 

a parametrized surface, i.e a smooth map F : ^ 

which is a homeomorphism onto its image and whose differential is injective. 

fu and fv 

derivatives of / with respect to the coordinates u and v respectively. 

T 

the set of parametrized surfaces in or pre-shape space. 

L{F{t)) 

length of the path t F{t) of parametrized surfaces. 

dist{Fi , F 2 ) 

infimum of the lengths of paths in F between Fi and F 2 , 

TfX 

the tangent space to the pre-shape space at a given parametrized shape /, 

iXl,X 2 )) 

a Riemannian metric on the pre-shape space T, 

^t{F\,F2) 

a geodesic path in F, from Fi to F 2 , parameterized by t G [0,1], 

00 (^1 ,F 2 ) = Fi, 01 (Fi, F 2 ) = F 2 , 

B/(Xi,X2) 

a symmetric bilinear map on the tangent space TfF which is non-negative 
when apply to the same tangent vector. 

KerS/ 

Kernel of the symmetric bilinear map Bf defined as 

KerBf = {Xi G TfF such that Bf{Xi,X 2 ) = 0, VX 2 G TfF}, 

((^ 1 ,^ 2 )) 

a non-negative semi-definite inner product on F, i.e. a symmetric bilinear map Bf 
on each tangent space TfF which is non-negative when apply to the same tangent vector. 

SO(3) 

the group of rotations in 

SO(3) XI 

the group of rotations and translations in 

Diff(S2) 

the group of all diffeomorphisms of S^, i.e. reparametrization group. 

r = Diff+(§2) 

the group of orientation-preserving diffeomorphisms of 
i.e. group of reparametrizations preserving orientation. 

G = Diff+(S2) X SO(3) X R3 

shape-preserving group. 

TjG 

quotient space of the pre-shape space by the group of reparametrizations preserving orientation, 
is identified with the shape space S by the application which maps 
a parametrized surface / : ^ to its image. 

[/] 

an element of the quotient space F/G which is the orbit of a parametrized shape F 
under the reparametrization group E = Diff+(S^) and the group SO(3) x 
[/] = {(^/ + 0 7 for 0 G SO(3), G M ^,7 G F} 




